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1. Introduction 

It is by now well established that the properties of matter formed from strongly interacting 
elementary particles change drastically at high temperatures. Quarks and gluons are no longer 
confined to move inside hadrons but organize in a new form of strongly interacting matter, the 
so-called quark-gluon plasma (QGP). The transition from hadronic matter to the QGP as well as 
properties of the high temperature phase have been studied extensively in lattice calculations over 
recent years [p. Nonetheless, detailed quantitative information on the transition and the structure 
of the high temperature phase in the physical situation of two light and a heavier strange quark 
((2 + l)-flavor QCD) is rare [^, ||, ||, §|. In order to relate experimental observables determined in 
relativistic heavy ion collisions to lattice results, it is important to achieve good quantitative control, 
in calculations with physical quark masses, over basic parameters that characterize the transition 
from the low to the high temperature phase of QCD. The most fundamental quantities characterize 
bulk properties of hot and dense matter: the transition temperature, energy density and pressure. 

Lattice QCD currently is the only quantitative approach to finite temperature QCD based on 
first principle calculation. At non-zero density however, lattice QCD is harmed by the sign problem 
ever since its inception. To overcome the sign problem is a challenging and outstanding problem. 
Nevertheless, during the last few years a lot of progress has been made to circumvent the sign 
problem for small values of pL q /T , where \i q is the quark chemical potential and T the temperature 

Mi 

2. Lattice formulation and calculational setup 

We study the thermodynamics of QCD with two light quarks (m/ = rh u = rhd) and a heavier 
strange quark (m s ) described by the QCD partition function which is discretized on a four dimen- 
sional lattice of size N% x N T , 

Z(p,m h m s ,N a ,Nj = f RdU^ (det D(m,)) 1/2 (det D(m s )) { ' A . (2.1) 

Here we will use staggered fermions to discretize the fermionic sector of QCD. The fermions have 
already been integrated out, which gives rise to the determinants of the staggered fermion matri- 
ces, D(rh[) and D(m s ) for the contributions of two light and one heavy quark degree of freedom, 
respectively. Moreover, j8 = 6/g 2 is the gauge coupling constant, m s j denote the dimensionless, 
bare quark masses in units of the lattice spacing a, and So is the gauge action which is expressed 
in terms of gauge field matrices U x ^ G St/ (3) located on the links (x,fl) = (xo,x,fi) of the four 
dimensional lattice; ;U = 0, ...,3. 

In our calculations we use a tree level, €?(a 2 ) improved gauge action, So, which includes the 
standard Wilson plaquette term and the planar 6-link Wilson loop. In the fermion sector, we use 
an improved staggered fermion action with 1-link and bended 3-link terms. The coefficient of the 
bended 3-link term has been fixed by demanding a rotationally invariant quark propagator up to 
0(p 4 ), which improves the quark dispersion relation at £?(a 2 ). This eliminates Gia 2 ) corrections 
to the pressure at tree level and leads to a strong reduction of cut-off effects in other bulk thermo- 
dynamic observables in the infinite temperature limit, as well as in G(g 2 ) perturbation theory [^]. 
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The 1-link term in the fermion action has been 'smeared' by adding a 3-link staple. This improves 
the flavor symmetry of the staggered fermion action Jlo|]. We call this action the p4fat3 action. It 
has been used previously in studies of QCD thermodynamics on lattices of temporal extent N t =4 
with larger quark masses [(2, 11]. 

Our studies of the transition to the high temperature phase of QCD [ |T2| ] have been performed 
on lattices of size N a X N T with N T = 4 and 6 and spatial lattice sizes N a = 8, 16, 24 and 32. We 
performed calculations for several values of the light to strange quark mass ratio, rhi/m s for fixed 
m s . The strange quark mass has been chosen such that the extrapolation to physical light quark 
mass values yields approximately the correct physical kaon mass value. The range of the light 
quark mass corresponds to a regime of the pseudo-scalar (pion) mass of 150 MeV</n ps <500 MeV 

In order to convert lattice units to physical units, zero temperature calculations are necessary 
which have been performed on 16 3 x 32 lattices. We use parameters characterizing the shape of 
the static quark potential (ro, r\, y/a) as well as hadron masses to set the scale for thermodynamic 
observables. 

The numerical simulation of the QCD partition function has been performed using the RHMC 
algorithm [13]. Unlike the hybrid-R algorithm [14] used in most previous studies of QCD thermo- 
dynamics performed with staggered fermions, this algorithm has the advantage of being exact, i.e. 
finite step size errors arising from the discretization of the molecular dynamics evolution of gauge 
fields in configuration space are eliminated through an additional Monte Carlo accept/reject step. 
This is possible with the introduction of a rational function approximation for roots of fermion 



determinants appearing in Eq. 2. 1 



3. Order parameters and susceptibilities 

To determine the QCD transition temperature and phase diagram, order parameter of the QCD 
transition are indispensable. In the chiral limit the chiral condensate (xj/yf) is the order parameter 
for the spontaneous chiral symmetry breaking of QCD. On the other hand in the heavy quark limit 
the Polyakov loop (L) is the order parameter of the deconfinement phase transition. For finite quark 
masses, these observables remain good indicators for the (pseudo) critical point. Especially their 
susceptibilities are useful to determine the pseudo critical coupling /3 C in numerical simulations. 

In Figure [T](left) we compare results for the light quark chiral condensate calculated on lat- 
tices of size 8 3 x 4 and 16 3 x 4. It clearly reflects the presence of finite volume effects at small 
values of the quark mass. While finite volume effects seem to be negligible for rhi/m s > 0.2, for 
rhi/m s = 0.1 we observe a small but statistically significant volume dependence for the chiral con- 
densate as well as for the Polyakov loop expectation value. This volume dependence is even more 
pronounced for rhi/m s = 0.05 and seems to be stronger at low temperatures. While the value of the 
chiral condensate increases with increasing volume the Polyakov loop expectation value decreases 
(Figure fright)). 

We use the Polyakov loop susceptibility as well as the disconnected part of the chiral suscep- 
tibility to locate the transition temperature to the high temperature phase of QCD, 

XL ^Nl((L 2 )-(L) 2 ) , (3.1) 
ff = ^3 (((TrD- 1 ^)) 2 )-^^- 1 ^)) 2 ) ,* = /,,. (3.2) 
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Figure 1: The light quark chiral condensate in units of a~ 3 (left) and the Polyakov loop expectation value 
(right) as function of the bare light quark mass in units of the temperature, mi/T = rhiNx for fixed /3 and 
m s = 0.065 on lattices of size 8 3 x 4 (circle) and 16 3 x 4 (triangles). Shown are results for various values of /3 
ranging from = 3.28 to j3 =3.4 (top to bottom for (yip) and bottom to top for (L)). Full and open symbols 
show results obtained from direct simulations and Ferrenberg-Swendsen interpolations, respectively. 
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Figure 2: The disconnected part of the light quark chiral susceptibility (left) and the Polyakov loop (right) 
on lattices of size 8 3 x 4 (squares) and 16 3 x 4 (circles) for two different values of the light quark mass. 
The curves show Ferrenberg-Swendsen interpolations of the data points obtained from multi-parameter his- 
tograms with an error band coming from Ferrenberg-Swendsen reweightings performed on different jack- 
knife samples. 

In Figure § we show results for the disconnected part of the light quark chiral susceptibility, 
Xi and the Polyakov loop, calculated on 8 3 x 4 and 16 3 x 4 lattices. The location of peaks in the 
susceptibilities has been determined from a Ferrenberg-Swendsen reweighting of data in the vicin- 
ity of the peaks. Errors on the critical couplings determined in this way have been obtained from 
a jackknife analysis where Ferrenberg-Swendsen interpolations have been performed on different 
sub-samples. In agreement with earlier calculations we find that the position of peaks in Xi an d 
Xl show only little volume dependence and that the peak height changes only little, although the 
maxima become somewhat more pronounced on the larger lattices. This is consistent with the 
transition being a crossover rather than a true phase transition in the infinite volume limit. 

Although differences in the critical coupling extracted from Xl and Xt are small we find that on 
small lattices the peak in the Polyakov loop susceptibility is located at a systematically larger value 
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of the gauge coupling /3 . In a finite volume this is, of course, not unexpected, and in the infinite 
volume limit an ambiguity in identifying the transition point may also remain for a crossover tran- 
sition. Nonetheless, we observe that the difference P C ,L — Pc,l decreases with increasing volume and 
is within errors consistent with zero for 16 3 x 4, which has the largest spatial volume expressed in 
units of the temperature, TV 1/3 = 4. On the smallest lattice, 8 3 x 4, we find j8 c ,l - (5 c j - 0.0077(9). 
Within the statistical accuracy of our data we also do not find any systematic quark mass depen- 
dence of this difference, j8 f ^ — j3 Cj ;. 

In addition to the light quark condensate and its susceptibility we also have analyzed the 
strange quark condensate and its susceptibility, % s . We find that the light and heavy quark con- 
densates are strongly correlated, which is easily seen in the MD-time evolution of these quantities. 
Already on the smallest lattices the position of the peak in the heavy quark susceptibilities is con- 
sistent with that deduced from the light quark condensate. On the larger, N a = 16, lattices the 
difference |j3 Ci / — j3 Cji | is in all cases zero within statistical errors, which are about 3 • 10~ 3 . Any 
temperature difference in the crossover behavior for the light and strange quark sector of QCD, 
which sometimes is discussed in phenomenological models, thus is below the 1 MeV level. 

We also have analyzed derivatives of the QCD partition function with respect to a quark chem- 
ical potential (see Section [7]). We note here that also from the analysis of the Taylor expansion 
coefficient C4 (J4) we find values for the pseudo-critical couplings that are in agreement with the 
above estimates [15] (see the quartic strange quark number fluctuations in Fig. [7] (left)). How- 
ever, as mentioned above the consistence of pseudo-critical couplings as determined from different 
observables is not necessary. If the transition becomes smoother closer to the continuum limit dif- 
ferent observables may lead to different estimates for the crossover point. In fact, for the stout-link 
improved staggered action and lattices with temporal extent of A^ T = 8 to 10, a large difference in f3 c 
determined from the chiral condensate and the strange quark number fluctuations has been found 



4. Scale setting and the heavy quark potential 



In order to calculate the transition temperature in terms of an observable that is experimentally 
accessible and can be used to set the scale for T c we have to perform a zero temperature calculation 
at the critical couplings f$ c determined in the previous section. This will allow us to eliminate 
the unknown lattice cut-off, a(j8 c ), which determines T c on a lattice with temporal extent N t , i.e. 
T c = l/N T a(P c ). To do so we have performed calculations at zero temperature, i.e. on lattices of 
size 16 3 x 32, and calculated several hadron masses as well as the static quark potential. From the 
latter we determine the string tension and extract short distance scale parameters rn, r\, which are 
defined as separations between the static quark anti-quark sources at which the force between them 



attains certain values [17], 



2<^(r) 



dr 



1.65 



r=r 



dr 



1.0 



(4.1) 



Although these scale parameters are not directly accessible to experiment they can be well esti- 
mated from heavy quarkonium phenomenology. Moreover, they have been determined quite accu- 



rately in lattice calculations through a combined analysis of the static quark potential [18] and level 



5 



QCD thermodynamics at zero and non-zero density 



Christian Schmidt 




Figure 3: Fit of the static quark potential with the Cornell ansatz and improved distances (left) and fit of the 
scale parameter ro versus the lattice coupling j3 to a renormalization group inspired ansatz (Eq 4.3) (right). 

splittings in bottomonium spectra [|l^]. Both these calculations have been performed on identical 
sets of gauge field configurations. We will use the value for ro determined in the bottomonium 
calculation [|^] for all conversions of lattice results to physical units, 

r = 0.469(7) fm . (4.2) 

Our zero temperature calculations have been performed at values of the gauge coupling in the 
vicinity of the j8 c we found from our analysis of the chiral susceptibilities. We typically generated 
several thousand configurations and analyzed the hadron spectrum and static quark potential on 
every 10 th configuration. We obtain the scale by using the simple Cornell form to fit our numerical 
results for the static quark potential, Vq q (r) = —a/r+or + c. With this fit- ansatz, which does 
not include a possible running of the coupling a, the force entering the definition of ro is easily 
calculated and we find from Eq. |4.l| , ro = \J (1.65 — a)jo . When fitting the potential, We replace 
the Euclidean distance on the lattice by an improved distance rj/a which relates the separation 
between the static quark and anti-quark sources to the Fourier transform of the tree-level lattice 
gluon propagator. This procedure removes most of the short distance lattice artifacts. We show an 
example for a fit of the static quark potential with improved distances in Fig. [3] (left). 

We have determined the scale parameter ro in units of the lattice spacing for 9 different param- 
eter sets (mi,m s ,P). This allows to interpolate between different values of the gauge coupling and 
quark masses. We use a renormalization group inspired ansatz [|(]] which takes into account the 
quark mass dependence of r^/a [J3j] and which approaches, in the weak coupling limit, the 2-loop 
j8-function for three massless flavors, 

{r /a)- 1 =R(p)(l +Ba 2 (p) + Ca\p))& A{m ' + ^ )+D . (4.3) 

Here /?(/3) denotes the 2-loop j3-function and a(j3) = R(p)/R(p) with j8 = 3.4 chosen as an ar- 
bitrary normalization point. We use this interpolation formula to set the scale for the transition 
temperature, a fit is shown in Fig. |3| (right). 
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Figure 4: T c ro (left) and T c /^fo (right) as a function of m ps ro on lattices with temporal extent N T = 4, 
m s = 0.065 (squares) and m s = 0.1 (triangles) as well as for N T = 6, m s = 0.04 (circles). Thin error bars 
represent the statistical and systematic error on ro/a and \faa. The broad error bar combines this error with 
the error on /3 C . The vertical line shows the location of the physical value m ps ro = 0.321(5) and its width 
represents the error on ro. The three parallel lines show results of fits based on Eq. 5.1 with d = 1.08 for 
N T = 4, 6 and Nz — ► °° (top to bottom). 

5. The transition temperature in (2+l)-flavor QCD 

To obtain the transition temperature we use the results for the scales ro/a and \[ca obtained 
from fits to the static quark potential. In cases where zero temperature calculations have not been 
performed directly at the critical coupling but at a nearby j8 -value we use Eq. |43| to determine 
the scales at p c (m/,m s ,N T ). The transition temperature is then obtained as T c tq = (ro/a)/N T or 
T c / \/a = 1 /y/aaNx- We show these results as function of the pseudo-scalar (pion) mass expressed 
in units of ro in Figure [| There we give 2 errors on T c ro and T c /y/o. A thin error bar reflects 
the combined statistical and systematic errors on the scales ro / a and \/aa obtained from our 3- 
parameter fit to the static quark potential. The broad error bar combines this uncertainty of the zero 
temperature scale determination with the scale-uncertainty arising from the error on j8 c . As can be 
seen, the former error, which typically is of the order of 2%, dominates our uncertainty on T c ro and 
T c j 'y/a on the coarser N t = 4 lattices, while the uncertainty in the determination of f5 c becomes 
more relevant for N t = 6. 

We have extrapolated our numerical results for T c tq and T c / yfo, which have been obtained for 
a specific set of lattice parameters (rhi,m s ,N x ), to the chiral and continuum limit using an ansatz that 
takes into account the quadratic cut-off dependence, (aT) 2 = 1 /N 2 , and a quark mass dependence 
expressed in terms of the pseudo-scalar meson mass, 

Y A ,AM = !W- +A(m ps r ) d + B/N 2 , Y = T c r 0l T c /^fc , (5.1) 

If the QCD transition is second order in the chiral limit the transition temperature is expected 
to depend on the quark mass as m) , or correspondingly on the pseudo-scalar meson mass as 
m 2 J/ S with d = 2/(58 ~ 1.08 characterizing universal scaling behavior in the vicinity of second 
order phase transitions belonging to the universality class of 0(4) symmetric, 3-dimensional spin 
models. If, however, the transition becomes first order for small quark masses, which is not ruled 
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out for physical values of the strange quark mass, the transition temperature will depend linearly 
on the quark mass (d = 2). A fit to our data set with d as a free fit parameter would actually favor 
a value smaller than unity, although the error on d is large in this case, d = 0.6(7). 

Fortunately, the extrapolation to the physical point is not very sensitive to the choice of d as 
our calculations have been performed close to this point. It does, however, increase the uncertainty 
on the extrapolation to the chiral limit. We have performed extrapolations to the chiral limit with d 
varying between d = 1 and d = 2. From this we find 

m ps r = 0: T c r = 0.444(6)+' 2 , T c /^ = 0.398(6)+}° , (5.2) 

where the central value is given for fits with the 0(4) exponent d = 1.08 and the lower and upper 
systematic error correspond to d = 1 and d = 2, respectively. Using the fit values for the parameter 
A that controls the quark mass dependence of T c r (A = 0.041(5)) and T c /^/a (A = 0.029(4)), 
respectively, we can determine the transition temperature at the physical point, fixed by m ps ror, 
where we then obtain a slightly larger value with reduced systematic errors, 

m ps r = 0.321(5) : T c r = 0.457(7)+| , T c /y/a = 0.408(8)1? . (5.3) 

Here the error includes the uncertainty in the value for the physical point, m ps ro, arising from the 
uncertainty in the scale parameter ro = 0.469(7) fm. We note that the extrapolated values for T c ro 
and T c / ^/a may also be interpreted as a continuum extrapolation of the shape parameters of the 
static potential. This yields r^yfc ~ 1.11 which is consistent with the continuum extrapolation 
obtained with the asqtad-action [J3]] . 



The fit parameter B which controls the size of the cut-off dependent term in Eq. qT| is in all 
cases close to 1 /3. We find B = 0.34(9) for fits to T c r and B = 0.33(7) for fits to T c /y/a, respec- 
tively. The critical temperatures for N x =4 thus are about 5% larger than the extrapolated value, 
and for N r = 6 the difference is about 2%. We therefore expect that any remaining uncertainties 
in our extrapolation to the continuum limit which may arise from higher order corrections in the 
cut-off dependence of T c tq are not larger than 2%. 

The results for the transition temperature obtained here for smaller quark masses and smaller 
lattice spacings is entirely consistent with the results for 2-flavor QCD obtained previously with the 



p4fat3 action on N z = 4 lattices in the chiral limit, T c /y/o = 0.425 ( 15 ) [JTT|] . We now find for (2+ 1 )- 
flavor QCD for N t = 4 in the chiral limit T c f \fo = 0.419(6). The continuum extrapolated result is, 
however, somewhat larger than the continuum extrapolated result obtained with the asqtad-action 
for (2+l)-flavor QCD in the chiral limit 1 , T c ro = 0.402(29) [Sj, which is based on the determination 
of transition temperatures on lattices with temporal extent N t = 4, 6 and 8. 

Although we frequently have referred to the physical value of ro during the discussion in the 
previous chapters we stress that our final result for dimensionless quantities, in particular T c vq and 



T c / \/o given in Eq. 5.2, does not depend on the actual physical value of ro or y/o. 



At present the scale parameter ro, deduced from the bottomonium level splitting using NRQCD 
[19], seems to be the best controlled lattice observable that can be used to set the scale for T c . Using 



'in [^] T c is given in units of r\ using results for r\ja taken from We have expressed T c in units of ro using 
= 1-4795 to convert r\ to the ro scale used by us. 
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for ro the value given in Eq. 4.2 we obtain for the transition temperature in QCD at the physical 
point, 

T c = 192(7) (4) MeV, (5.4) 



where the statistical error includes the errors given in Eq. 5.3 as well as the uncertainty in the value 
of ro an d the second error reflects our estimate of a remaining systematic error on the extrapola- 



tion to the continuum limit. As discussed after Eq. we estimate this error which arises from 



neglecting higher order cut-off effects in our ansatz for the continuum extrapolation, Eq. [5JJ to be 
about 2%. Our resut for T c agrees within our current uncertainties with the result by the MILC 
Collaboration which is T c = 169(12) (4) [S]. However, we encounter a discrepancy with the result 
by the Wuppertal group of T c = 15 1(3) (3) [ |I6[ ] as determined by the chiral condensate. 



6. The equation of state 

We calculate the equation of state on a line of constant physics (LCP). Along this line the phys- 
ical quark masses are fixed but the bare light and strange quark masses, which are the parameters in 
our action, change as a function of the coupling. We define the LCP by demanding that the ratio of 
strange pseudo-scalar mass over kaon mass (m^/m^) as well as m Ss in units of the scale parameter 
(roniss) stay constant. We find that in good approximation the first condition is fulfilled by holding 
the ratio of bare light and strange quark mass (m//m v ) fixed. We calculate meson masses in a wide 
range of the parameter space and use a renormalization group inspired interpolation formula to 
carefully tune the light quark masses as a function of the coupling /3 in order to satisfy the second 
condition. 

For calculating the pressure and the interaction measure we employ the integral method. Since 
the logarithm of the partition function can not be calculated easily on the lattice, derivatives with 
respect to the bare parameters are calculated (dlnZ/dp, dlnZ/dmi, d\nZ/dm s ) and integrated 
along the LCP. We find for the pressure 



P_ 
y4 



ft 



K dp 1 

JPo 



1 



N*N T 



((S g )q-(Sg)t) 



2((w)/o 
mi({WW) s o-(YY) sT 



(W)l T ) + — (( W)sO - (W) sT ] 

mi 



dp' 



th s jyh\ 



drhjriii 
dp' 



(6.1) 



As we find that rii s /riii = 10 can be kept fixed on a LCP, the last term in the integral vanishes. To 
calculate the interaction measure e — 3p, the lattice p -function has to be known. Here we have 

dp/T* 



3p 



dr Vr 4 

e-3p 

y4 



da 



+ 



gluon 



LCP 

3p y 



dp 



rp4 



fermion 



3p 



(6.2) 



m s jmi 



Again the last term vanishes to good accuracy on the LCP. We emphasize that a precise knowledge 
of the P -function, adp/da, along the LCP is necessary to calculate the interaction measure and 
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Figure 5: The interaction measure (left) and pressure (right) along a line of constant physics for N z — 4 and 
6. The ratio of light to strange quark mass has been held fixed at m//m s = 0. 1 . Also shown are earlier results 
for the pressure, obtained with 3-flavor QCD and const. m q /T = 0.4 [2]. 

energy density. We calculate the /3 -function from measurements of the scale parameter vq, which 
we perform on each temperature (coupling) which we include in our equation of state. 

In Fig. H we show our preliminary results for the interaction measure (left) and the pressure 
(right) on the LCP for N t =4 and 6. The ratio of light to strange quark mass has been held fixed 
at rhi/m s = 0.1. We compare our result for the pressure with earlier results for 3-flavor QCD and 
a constant m q /T = 0.4 [0]. One clearly sees a mass dependence, however, the mass dependence is 
small. Note that the quark masses differ by more than an order of magnitude in the high temperature 
region. 

Or results on the interaction measure and pressure obtained with the p4fat3 improved action 



are in complete agreement with corresponding results obtained with asqtad fermions [|22|]. Further- 
more, the results show little cut-off dependence in the entire temperature regime analyzed. This 
is in contrast to calculations performed with the standard staggered discretization scheme [|Z3]], 
which leads to large cut-off effects in the hight temperature limit for lattices with temporal extent 
of N T = 4 and 6. 

7. Hadronic fluctuations at zero and non-zero density 

It is conceptually very simple to calculate the expansion coefficients of any observable O in a 
Taylor series around pL q = where \i q is the quark chemical potential: 

O{(L)=c + c l (L+ l -c 2 (L 1 + --- . (7.1) 

Since on the lattice all quantities are given in units of the lattice spacing (a), the expansion param- 
eter is jx = a\l q = N~ l (jjL q /T). This idea goes back to the first calculation of the quark number 



susceptibility [|24|]. The response of hadron masses [25] as well as the pressure and further bulk 



thermodynamic quantities [[7],|2(], 27, ^8|] have been studied by this method. The first two nontrivial 



coefficients in Eq. (7.1) are given by 



,'dO\ I dlndetD\ 
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/d 2 0\ /d<9dlndetD\ / d 2 lndetD\ / \ /d 2 lndetD\ 

Besides derivatives of the observable itself, the calculation of derivatives of IndetD with respect to 
jl is required. The derivatives have to be taken at jX = 0. Note that due to the particle-antiparticle 
symmetry of the partition function (Z(n q ) = Z(—jj, q )) all odd coefficients in Eq. ( fP| ) vanish iden- 
tically. For the same reason we have (dlndetD/d/t) = at jx = 0. We explicitly use this property 



in Eq. ( |7.2| ) to derive the expansion coefficients. 

The advantages of this method are that expectations values only have to be evaluated at jx = 0, 
i.e. calculations are not directly affected by the sign problem. Furthermore, all derivatives of the 
fermion determinant can be expressed in terms of traces by using the identity IndetD = TrlnD. 
This enables the stochastic calculation of the expansion coefficients by the random noise method, 
which is much faster than a direct evaluation of the determinant. Moreover, the continuum and 
infinite volume extrapolations are well defined on a coefficient by coefficient basis. 

Quark number fluctuations % q belong to the most important observables that allow to follow 
the transition line into the non-zero chemical potential plane. They diverge at the critical end-point 
and thus provide an excellent signal for its existence and location. Eventually these fluctuations 
may be detectable in heavy ion experiments. Hadronic fluctuations can be computed from Taylor 
expansion coefficients of the pressure with respect to the quark chemical potential: 

P £ /„s / £i« \ " ., 1 AT? dlnZ 



with c n (T) 



(7.3) 



T 4 to KT J n\Nid(fLN T )« 

Note that in the Taylor expansion of the pressure the up and down quark chemical potentials have 
been chosen to be equal. Having calculated the coefficients c n {T) one can construct the quark 
number density and quark number fluctuations 

^=inCn{T)(f)~ 1 ; ff = £n(n-l)c„(r)(^)" 2 . (7.4) 

In the case of two flavors of p4-improved staggered fermions, with m q /T = 0.4 the first three 
non-zero coefficients c%, cn, and C(, have been calculated [27] and are shown in Fig. |6| We also 



show our preliminary results on C2 for (2+l)-flavor QCD and a quark mass ratio of m q /m s = 0.1. 
A mass dependence is clearly evident. 

In Fig. [7] (left) we show the quartic quark number fluctuations of strange quarks from our 
ongoing (2+l)-flavor QCD simulations, which are given by rff = (\/{VT 3 ))(d A \nZ/d(}i s /T) A ). 
As one can see also the quartic strange quark fluctuations show a peak at T c . 

In Fig. fright) we show the quark number fluctuations for tif = 2 and m q /T = 0.4 [27] for 



several values of the quark chemical potential including only the leading (mu q /T) 2 correction 
which is proportional to C4 . It is interesting to see that at \i q = 0, % q shows a rapid but monotonic 
increase at the transition temperature, whereas a cusp is developing at T c {jx q ) for ji q > 0. This is a 
clear sign for approaching the critical end-point. 

8. The critical end-point 

Locating the critical point is one of the most challenging goals of lattice QCD calculations at 
finite chemical potential. The first attempt to locate the critical point used the reweighting method 
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Figure 6: The Taylor expansion coefficients c%, C4 and eg of the pressure for fit = 2 and m q /T = 0.4 [27]. 
Also shown are preliminary results for (2+l)-fiavor QCD and a quark mass ratio of m q /m s =0.1. 




Figure 7: Quartic strange quark number fluctuations for (2+l)-flavor QCD and a quark mass ratio of 
ihq/m s = 0.1 (left) and the quark number susceptibility for nj — 2 and m q /T = 0.4 at several values of 
the quark chemical potential (right). 



[29]. For this calculation, 2+1 flavor of standard staggered fermions have been used at a pion 



mass of about 300 MeV and a kaon mass of about 500 MeV. Lattice sizes have, however, been 
rather small (4 3 x 4 - 8 3 x 4). A critical chemical potential of ]X% n = 725(35) MeV was found. A 
second calculation using again the reweighting method, with physical masses (m n = 150 MeV, 
m K = 500 MeV) and somewhat larger volume (6 3 x 4 - 12 3 x 4), let to }X c B rit = 360(4) MeV 

When using the reweighting method for locating the critical point, the minima of the normal- 
ized partition function in the complex j3 -plane (Lee-Yang zeros) have to be determined 



7 = 



Z(j8 Re ,j3i m ,Ai) 



Z(/3 Re ,0,0) 



6ij3iV T iv3AS Ge (e e (JV / /4)(liidetD( J u)-lndetD(0))\ 



(Aic ,0,0) 



(8.1) 



In SU (3) gauge theory, where we have e = 1, this can be done with high accuracy [30]. In QCD 



with non-zero chemical potential the analysis of Lee- Yang zeros is, however, subtle [30]. For large 
volumes and chemical potentials the phase factor of the determinant e' 6 will force the Lee- Yang 
zero onto the real axis, which might lead to an underestimation of the critical point. 
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Figure 8: Estimates of the radius of convergence in the (T,jj, q ) -plane (left), the ratio C4/C2 of the expansion 
coefficients (right). The horizontal lines indicate the resonance gas limit (T — > 0) and the SB limit in the 
continuum and at N T = 4 (T — > °°). 



Another difficulty with the reweighting method at finite chemical potential has been pointed 
out in [31]. It was noted, that taking the fourth (or square) root of the determinant (which is 
necessary in order to simulate 2 or 1-flavor QCD with staggered fermions; see also [f32|]) could lead 
to phase ambiguities. This problem becomes acute when \i q > m n /2. 

All of the above mentioned limitations are, however, irrelevant for the location of the critical 
point with the reweighting method if the critical point is located at small values of }i q . 

Using the Taylor expansion coefficients of the pressure, it is also possible to estimate the 
location of the critical point. The convergence radius of the expansion is limited by the nearest 
singularity in the complex chemical potential plane. For each fixed temperature, the radius of 
convergence is given by 



lim p n 




Cn+2 



(8.2) 



Moreover, the sign of the coefficients c n gives information about the location of the singularity 
in the complex plane. If all coefficients are positive, the singularity is located on the real axis of 
the complex chemical potential plane. If the sign is strictly alternating, the singularity lies on the 
imaginary axis. For a detailed discussion see [33]. 

Having only a limited number of expansion coefficients, one can only estimate p. The hope 
is that the convergence of the p„ will be fast. Indeed, a clustering of the p n is seen in the phase 
diagram, as shown in Fig. [8] [27]. This calculation, which has been performed with 2 flavors of p4 
improved fermions and m K /m p = 0.7, suggests a critical chemical potential of jig" 500 MeV. 
All calculated p„ are, however, consistent (within statistical error) with the resonance gas model in 
the Boltzmann approximation, where the radius of convergence is infinity. 

The authors of fl34] ] have estimated the critical chemical potential from a Taylor expansion 
of the quark number susceptibility and find /Xg" 180 MeV. Two flavors of standard staggered 
fermions have been used on lattices up to 24 3 x 4 and quark mass corresponding to m„/m p = 0.3. 
The difference between the two estimates [21 , 34 ] of the critical point is large. We note that the 
second estimate comes from the expansion coefficients of % q . As can be seen from Eq. 7.4 this 
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Figure 9: Phase diagram from the density of state method [35]. 



will result in a smaller p n for each fixed n. The limit lim„^oo p n is of course the same. For finite n, 
however, the estimate of pL^" rj 1 80 MeV would correspond to pL^" rj 240 MeV, when estimating 
the p n with coefficients of the same order from the expansion of the pressure. Nonetheless, the 
difference between the two estimates is still striking. The origin could be the difference in mass. 
However, preliminary results from the RBC-Bielefeld Collaboration, also shown in Fig. [| do not 
indicate a strong mass dependence in 04/02 = l/pf- 



9. Beyond the critical point 

Even more challenging than locating the critical point, is the study of the physics at high den- 
sities and low temperatures. One attempt to do so is a calculation using the density of states method 



[35]. Using four flavors of standard staggered fermions (i.e. taking the root of the determinant is 



not necessary), several simulation points in the (j3,/t) plane have been chosen to generate phase 



quenched configurations by employing the method proposed in [36]. The lattice size has been 
6 3 x 4, 6 4 and 6 3 x 8. The quark mass was chosen to be m/T = 0.3. The generation has been 
done with constrained plaquettes. In oder to do so, we introduce a sharply peaked Gaussian po- 
tential in to the partition function, which in practice leads to a modification of the force term of 
the HMD-R algorithm. For each simulation point, several runs have been performed with about 20 
different values of the plaquette. By calculating the eigenvalues of the reduced matrix the phase of 
the determinant was calculated for each of those runs. By numerically calculating the integrals 

(p) = j dxxp(x){cos(d)) x (P 2 ) = J dxx 2 p(x) {cos(d)) x , (9.1) 

we recover the grand canonical expectation value of the plaquette and its square. Here p (x) is the 
density of states, which has been measured by the integral method, usually used to calculate the 
pressure. The susceptibility of the plaquette is then given by the usual expression %p = (P 2 ) — (P) 2 - 
From the peak position of the plaquette susceptibility the phase diagram was calculated as shown 
in Fig. |9| The points at T = 93 MeV are from calculations on 6 3 x 8 lattices. Note, that we make 
no statement about the order of the transition lines. To determine the order of the transition one has 
to perform a finite-size-scaling analysis which is beyond the scope of this article. 
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The plaquette expectation value and plaquette susceptibility suggest three different phases, 
which coincide in a triple point. The triple point is located around pL q l rj 300 MeV, however its 
temperature (T tri ) decreases from T tri m 148 MeV on the 4 4 lattice to T tri rj 137 MeV on the 6 4 
lattice. This shift reflects the relatively large cut-off effects one faces, with standard staggered 
fermions and temporal extents of 4 and 6. 

The new phase at large chemical potentials and low temperatures is a natural candidate for a 
color superconducting phase. Recently, by combining experimental results from cold atoms in a 
trap j37| ] and some universal arguments, an upper bound for the transition line from the quark gluon 
plasma phase (QGP) to the superconducting phase (SC) was proposed (T c < 0.35E F ) [@|3|]. To 
first approximation the Fermi -Energy Ep is given by the chemical potential \i q . In p9^ ] the triple 
point was estimated by comparing this upper bound with the experimental freeze-out curve. A 
temperature of T tn < 70 MeV was found. Our value of the triple-point roughtly corresponds to 
T c < 0A6Ep. It is interesting that the two values are close. 

At low temperatures we find a phase boundary which is very steep and almost independent of 
jX q . Although our lowest temperature is 96 MeV an extrapolation to T = seems to be reasonable 
and would yield a critical chemical potential of pL q (T = 0) « 250 MeV or equivalently Hb/T c (ij.b = 
0) Ri 4.7. This number appears to be at the lower edge of the phenomenological expectation of 
IXb/T c {}Xb = 0) ~ 5 — 10. Note, that our lattice spacing is close to the strong coupling regime and 
we should feel the influence of the strong coupling limit. Strong coupling expansion calculations 
in general yield much lower values of }JLb/T c (pLb = 0)<1.5 p0|]. 

For this work the density of state method has been employed, which works well on small 
lattices up to chemical potentials of [J. q /T<3 (other methods [§|,[7[,[8|] worked up to ll q /T<Y). The 
method is however extremely expensive and thus will in the near future not yield results close to 
the thermodynamic limit or the continuum limit, due to limitations in computer resources. 

We have to emphasize that this simulations have been carried out on coarse lattices with an 
unphysical value of rif ■ = 4 degenerate fermion flavor, and that neither the continuum nor the ther- 
modynamic limit has been taken. Since we used unimproved staggered fermions, the corrections 
due to a finite lattice spacing are large. We also expect corrections due to the finite size of our 
volume. The simulations have not been performed with a constant quark mass, but m q /T = 0.3 has 
been held fixed. 
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